clear all
global dir = "C:\Users\mrueda\Documents\Emory\Papers\Networks_persistance\do_files\do_files_APSA21\post_JOP\Replication_BJPS\"	



	global outcomes2 fcontract fgot_above_ext fruns_any
	global outcomes3 donate_15any b5
	global outcomes4 fdonate_15any fb5 
	global outcomes5 nfdonate_15any nfb5

	
	*Load data

	
	********************************************
	*	DISCONTINUITY FIGURES
	********************************************

	cd "$dir\Figures"

	*Define the level of confidence
	local l = 95
	
	foreach x in nfcontract {
		use  "$dir\Data\cand_level_persist_rep.dta",clear
		drop if `x'==.
		// Basic estimation and plot with degreee 1
		rdrobust `x' margin_victory, vce(cluster muni_code) p(1) bwselect(mserd) level(`l') 
		local bw=round(`e(h_l)',0.001)
		
		// Gen weights
		gen weights = .
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory < 0 & margin_victory >= -`bw'
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory >= 0 & margin_victory <= `bw'

		rdplot `x' margin_victory, vce(cluster muni_code) p(1) nbins(24 24) genvars h(`bw' `bw') scale(2 2) kernel(triangular)
		// Not very nice with local linear
		gen bin_no=rdplot_id
		gen meanx_bin=rdplot_mean_x
		gen mean_ybin=rdplot_mean_y
		sort bin_no margin_victory
		bys bin_no: gen pos=_n
		replace mean_ybin=. if pos!=1  // Keep the first value per bin to avoid pltting several equal means
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) (lpolyci `x' margin_victory if margin_victory<0 & abs(margin_victory)<=`bw', level(`l') kernel(triangle) bw(`bw') deg(1) fcolor(none)) (lpolyci `x' margin_victory if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') bw(`bw') kernel(triangle) deg(1) fcolor(none)), xline(0)  legend(off)
		
		// Global linear as rdplot does
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) /// 
					 (lfitci `x' margin_victory [aw = weights] if margin_victory<0 & abs(margin_victory)<=`bw', ciplot(rline) level(`l') fcolor(none))  ///
					 (lfitci `x' margin_victory [aw = weights] if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') ciplot(rline) fcolor(none)), ///
					 xline(0)  legend(off)   graphregion(fcolor(white)) ///
					 xtitle(Margin Victory) 
					 
		graph export "Fig_D3_`x'.pdf", replace
		drop rdplot_id - pos weights
	}

		foreach x in $outcomes2 {
		use  "$dir\Data\cand_level_persist_rep.dta",clear
		drop if `x'==.
		// Basic estimation and plot with degreee 1
		rdrobust `x' margin_victory, vce(cluster muni_code) p(1) bwselect(mserd) level(`l') 
		local bw=round(`e(h_l)',0.001)
		
		// Gen weights
		gen weights = .
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory < 0 & margin_victory >= -`bw'
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory >= 0 & margin_victory <= `bw'

		rdplot `x' margin_victory, vce(cluster muni_code) p(1) nbins(24 24) genvars h(`bw' `bw') scale(2 2) kernel(triangular)
		// Not very nice with local linear
		gen bin_no=rdplot_id
		gen meanx_bin=rdplot_mean_x
		gen mean_ybin=rdplot_mean_y
		sort bin_no margin_victory
		bys bin_no: gen pos=_n
		replace mean_ybin=. if pos!=1  // Keep the first value per bin to avoid pltting several equal means
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) (lpolyci `x' margin_victory if margin_victory<0 & abs(margin_victory)<=`bw', level(`l') kernel(triangle) bw(`bw') deg(1) fcolor(none)) (lpolyci `x' margin_victory if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') bw(`bw') kernel(triangle) deg(1) fcolor(none)), xline(0)  legend(off)
		
		// Global linear as rdplot does
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) /// 
					 (lfitci `x' margin_victory [aw = weights] if margin_victory<0 & abs(margin_victory)<=`bw', ciplot(rline) level(`l') fcolor(none))  ///
					 (lfitci `x' margin_victory [aw = weights] if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') ciplot(rline) fcolor(none)), ///
					 xline(0)  legend(off)   graphregion(fcolor(white)) ///
					 xtitle(Margin Victory) 
					 
		graph export "Fig_D3_`x'.pdf", replace
		drop rdplot_id - pos weights
	}
	

		

	*Define the level of confidence
	local l = 95
	
	foreach x in $outcomes3 {
		use  "$dir\Data\cand_level_persist_rep.dta",clear
		drop if `x'==.
		// Basic estimation and plot with degreee 1
		rdrobust `x' margin_victory, vce(cluster muni_code) p(1) bwselect(mserd) level(`l') 
		local bw=round(`e(h_l)',0.001)
		
		// Gen weights
		gen weights = .
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory < 0 & margin_victory >= -`bw'
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory >= 0 & margin_victory <= `bw'

		rdplot `x' margin_victory, vce(cluster muni_code) p(1) nbins(24 24) genvars h(`bw' `bw') scale(2 2) kernel(triangular)
		// Not very nice with local linear
		gen bin_no=rdplot_id
		gen meanx_bin=rdplot_mean_x
		gen mean_ybin=rdplot_mean_y
		sort bin_no margin_victory
		bys bin_no: gen pos=_n
		replace mean_ybin=. if pos!=1  // Keep the first value per bin to avoid pltting several equal means
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) (lpolyci `x' margin_victory if margin_victory<0 & abs(margin_victory)<=`bw', level(`l') kernel(triangle) bw(`bw') deg(1) fcolor(none)) (lpolyci `x' margin_victory if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') bw(`bw') kernel(triangle) deg(1) fcolor(none)), xline(0)  legend(off)
		
		// Global linear as rdplot does
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) /// 
					 (lfitci `x' margin_victory [aw = weights] if margin_victory<0 & abs(margin_victory)<=`bw', ciplot(rline) level(`l') fcolor(none))  ///
					 (lfitci `x' margin_victory [aw = weights] if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') ciplot(rline) fcolor(none)), ///
					 xline(0)  legend(off)   graphregion(fcolor(white)) ///
					 xtitle(Margin Victory) 
					 
		graph export "Fig_D1_`x'.pdf", replace
		drop rdplot_id - pos weights
	}
	



	*Define the level of confidence
	local l = 95
	
		
	foreach x in $outcomes5 {
		use  "$dir\Data\cand_level_persist_rep.dta",clear
		drop if `x'==.
		// Basic estimation and plot with degreee 1
		rdrobust `x' margin_victory, vce(cluster muni_code) p(1) bwselect(mserd) level(`l') 
		local bw=round(`e(h_l)',0.001)
		
		// Gen weights
		gen weights = .
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory < 0 & margin_victory >= -`bw'
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory >= 0 & margin_victory <= `bw'

		rdplot `x' margin_victory, vce(cluster muni_code) p(1) nbins(24 24) genvars h(`bw' `bw') scale(2 2) kernel(triangular)
		// Not very nice with local linear
		gen bin_no=rdplot_id
		gen meanx_bin=rdplot_mean_x
		gen mean_ybin=rdplot_mean_y
		sort bin_no margin_victory
		bys bin_no: gen pos=_n
		replace mean_ybin=. if pos!=1  // Keep the first value per bin to avoid pltting several equal means
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) (lpolyci `x' margin_victory if margin_victory<0 & abs(margin_victory)<=`bw', level(`l') kernel(triangle) bw(`bw') deg(1) fcolor(none)) (lpolyci `x' margin_victory if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') bw(`bw') kernel(triangle) deg(1) fcolor(none)), xline(0)  legend(off)
		
		// Global linear as rdplot does
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) /// 
					 (lfitci `x' margin_victory [aw = weights] if margin_victory<0 & abs(margin_victory)<=`bw', ciplot(rline) level(`l') fcolor(none))  ///
					 (lfitci `x' margin_victory [aw = weights] if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') ciplot(rline) fcolor(none)), ///
					 xline(0)  legend(off)   graphregion(fcolor(white)) ///
					 xtitle(Margin Victory) 
					 
		graph export "Fig_D5_`x'.pdf", replace
		drop rdplot_id - pos weights
	}
	

		
	foreach x in $outcomes4 {
		use  "$dir\Data\cand_level_persist_rep.dta",clear
		drop if `x'==.
		// Basic estimation and plot with degreee 1
		rdrobust `x' margin_victory, vce(cluster muni_code) p(1) bwselect(mserd) level(`l') 
		local bw=round(`e(h_l)',0.001)
		
		// Gen weights
		gen weights = .
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory < 0 & margin_victory >= -`bw'
		replace weights = (1 - abs(margin_victory / `bw')) if margin_victory >= 0 & margin_victory <= `bw'

		rdplot `x' margin_victory, vce(cluster muni_code) p(1) nbins(24 24) genvars h(`bw' `bw') scale(2 2) kernel(triangular)
		// Not very nice with local linear
		gen bin_no=rdplot_id
		gen meanx_bin=rdplot_mean_x
		gen mean_ybin=rdplot_mean_y
		sort bin_no margin_victory
		bys bin_no: gen pos=_n
		replace mean_ybin=. if pos!=1  // Keep the first value per bin to avoid pltting several equal means
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) (lpolyci `x' margin_victory if margin_victory<0 & abs(margin_victory)<=`bw', level(`l') kernel(triangle) bw(`bw') deg(1) fcolor(none)) (lpolyci `x' margin_victory if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') bw(`bw') kernel(triangle) deg(1) fcolor(none)), xline(0)  legend(off)
		
		// Global linear as rdplot does
		tw (scatter mean_ybin margin_victory if abs(margin_victory)<=`bw', mcolor(gs10)) /// 
					 (lfitci `x' margin_victory [aw = weights] if margin_victory<0 & abs(margin_victory)<=`bw', ciplot(rline) level(`l') fcolor(none))  ///
					 (lfitci `x' margin_victory [aw = weights] if margin_victory>=0 & abs(margin_victory)<=`bw', level(`l') ciplot(rline) fcolor(none)), ///
					 xline(0)  legend(off)   graphregion(fcolor(white)) ///
					 xtitle(Margin Victory) 
					 
		graph export "Fig_D5_`x'.pdf", replace
		drop rdplot_id - pos weights
	}
	
	********************************************
	*	BANDWIDTH FIGURES
	********************************************

	cd "$dir"
	use "Data\cand_level_persist_rep.dta",clear
	cd "$dir\Figures"
	**Here is where the running variable goes.
	foreach x in margin_victory {
	**Here is where the outcome goes.
	*foreach y in $outcomes {
		foreach y in nfcontract {
		************
		* Linear
		*************
		*mserd approach, linear (get optimal bw)
		rdrobust `y' `x', all vce(cluster muni_code) p(1) level(95) 

		local bw=round(`e(h_r)',0.001)
		local bw_double = `bw'*2
		local bw_half = `bw'/2
		local rho =  `e(h_r)'/`e(b_r)'
		local counter=1
		local b_`counter' = round(`e(tau_cl)',0.001)
		local uci_`counter' = round(`e(ci_r_rb)',0.001)
		local lci_`counter' = round(`e(ci_l_rb)',0.001)
		local bw_`counter' = round(`e(h_r)',0.001)


		*According to Cattaneo, Idrobo and Titiunik, 
		*we should do double of the optimal bandwidht.
		local j = 1
		while `j'>0 {
			local counter = `counter'+1
			
			*CTT approach
			local b = `bw_double'/`rho'
			rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')
			
			local b_`counter' = round(`e(tau_cl)',0.001)
			local uci_`counter' = round(`e(ci_r_rb)',0.001)
			local lci_`counter' = round(`e(ci_l_rb)',0.001)
			local bw_`counter' = round(`e(h_r)',0.001)
			local bw_double = `bw_double'- (0.005)
			
			if `bw_double'<=`bw_half' {
				local j = 0
				
				* CTT approach
				local counter = `counter'+1
				local b = `bw_double'/`rho'
				rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')

				local b_`counter' = round(`e(tau_cl)',0.001)
				local uci_`counter' = round(`e(ci_r_rb)',0.001)
				local lci_`counter' = round(`e(ci_l_rb)',0.001)
				local bw_`counter' = round(`e(h_r)',0.001)
					
			
			}
			else {
				local j = 1
			}   
	}


	matrix graph = J(`counter',4,.)
	forvalues i=1/`counter' {
		mat graph[`i',1] = `bw_`i''
		mat graph[`i',2] = `b_`i''
		mat graph[`i',3] = `lci_`i''
		mat graph[`i',4] = `uci_`i''
		
	}

	mata : st_matrix("graph", sort(st_matrix("graph"), 1))
	preserve
	clear 
	svmat graph
	* graphs
	local line_2 = round(`x_line2',0.001)
	*attacks
	graph tw (sc graph2 graph1) (rcap graph4 graph3 graph1, lcolor(dkgreen)), title("") subtitle("") ///
	xtitle("Margin Victory") ytitle("Point Estimate") xli(`bw') xlabel(#8) xscale(range(`bw_half' `bw_double')) legend(off)  graphregion(fcolor(white))  yline(0)
	
	*Save graph
	graph export "Fig_E2_`y'.pdf", as(pdf) replace
	
	restore
	}
	}
	cd "$dir"
	use "Data\cand_level_persist_rep.dta",clear
	cd "$dir\Figures"
	foreach y in $outcomes2 {
		************
		* Linear
		*************
		*mserd approach, linear (get optimal bw)
		rdrobust `y' margin_victory, all vce(cluster muni_code) p(1) level(95) 

		local bw=round(`e(h_r)',0.001)
		local bw_double = `bw'*2
		local bw_half = `bw'/2
		local rho =  `e(h_r)'/`e(b_r)'
		local counter=1
		local b_`counter' = round(`e(tau_cl)',0.001)
		local uci_`counter' = round(`e(ci_r_rb)',0.001)
		local lci_`counter' = round(`e(ci_l_rb)',0.001)
		local bw_`counter' = round(`e(h_r)',0.001)


		*According to Cattaneo, Idrobo and Titiunik, 
		*we should do double of the optimal bandwidht.
		local j = 1
		while `j'>0 {
			local counter = `counter'+1
			
			*CTT approach
			local b = `bw_double'/`rho'
			rdrobust `y' margin_victory, all vce(cluster muni_code) level(95) h(`bw_double') b(`b')
			
			local b_`counter' = round(`e(tau_cl)',0.001)
			local uci_`counter' = round(`e(ci_r_rb)',0.001)
			local lci_`counter' = round(`e(ci_l_rb)',0.001)
			local bw_`counter' = round(`e(h_r)',0.001)
			local bw_double = `bw_double'- (0.005)
			
			if `bw_double'<=`bw_half' {
				local j = 0
				
				* CTT approach
				local counter = `counter'+1
				local b = `bw_double'/`rho'
				rdrobust `y' margin_victory, all vce(cluster muni_code) level(95) h(`bw_double') b(`b')

				local b_`counter' = round(`e(tau_cl)',0.001)
				local uci_`counter' = round(`e(ci_r_rb)',0.001)
				local lci_`counter' = round(`e(ci_l_rb)',0.001)
				local bw_`counter' = round(`e(h_r)',0.001)
					
			
			}
			else {
				local j = 1
			}   
	}


	matrix graph = J(`counter',4,.)
	forvalues i=1/`counter' {
		mat graph[`i',1] = `bw_`i''
		mat graph[`i',2] = `b_`i''
		mat graph[`i',3] = `lci_`i''
		mat graph[`i',4] = `uci_`i''
		
	}

	mata : st_matrix("graph", sort(st_matrix("graph"), 1))
	preserve
	clear 
	svmat graph
	* graphs
	local line_2 = round(`margin_victory_line2',0.001)
	*attacks
	graph tw (sc graph2 graph1) (rcap graph4 graph3 graph1, lcolor(dkgreen)), title("") subtitle("") ///
	xtitle("Margin Victory") ytitle("Point Estimate") xli(`bw') xlabel(#8) xscale(range(`bw_half' `bw_double')) legend(off)  graphregion(fcolor(white))  yline(0)
	
	*Save graph
	graph export "Fig_E2_`y'.pdf", as(pdf) replace
	
	restore
	}
	
	cd "$dir"
	use "Data\cand_level_persist_rep.dta",clear

	cd "$dir\Figures"
	******************
	*	FIGURES TABLE 3
	******************
	**Here is where the running variable goes.
	foreach x in margin_victory {
	**Here is where the outcome goes.
	*foreach y in $outcomes {
		foreach y in $outcomes3 {
		************
		* Linear
		*************
		*mserd approach, linear (get optimal bw)
		rdrobust `y' `x', all vce(cluster muni_code) p(1) level(95) 

		local bw=round(`e(h_r)',0.001)
		local bw_double = `bw'*2
		local bw_half = `bw'/2
		local rho =  `e(h_r)'/`e(b_r)'
		local counter=1
		local b_`counter' = round(`e(tau_cl)',0.001)
		local uci_`counter' = round(`e(ci_r_rb)',0.001)
		local lci_`counter' = round(`e(ci_l_rb)',0.001)
		local bw_`counter' = round(`e(h_r)',0.001)


		*According to Cattaneo, Idrobo and Titiunik, 
		*we should do double of the optimal bandwidht.
		local j = 1
		while `j'>0 {
			local counter = `counter'+1
			
			*CTT approach
			local b = `bw_double'/`rho'
			rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')
			
			local b_`counter' = round(`e(tau_cl)',0.001)
			local uci_`counter' = round(`e(ci_r_rb)',0.001)
			local lci_`counter' = round(`e(ci_l_rb)',0.001)
			local bw_`counter' = round(`e(h_r)',0.001)
			local bw_double = `bw_double'- (0.005)
			
			if `bw_double'<=`bw_half' {
				local j = 0
				
				* CTT approach
				local counter = `counter'+1
				local b = `bw_double'/`rho'
				rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')

				local b_`counter' = round(`e(tau_cl)',0.001)
				local uci_`counter' = round(`e(ci_r_rb)',0.001)
				local lci_`counter' = round(`e(ci_l_rb)',0.001)
				local bw_`counter' = round(`e(h_r)',0.001)
					
			
			}
			else {
				local j = 1
			}   
	}


	matrix graph = J(`counter',4,.)
	forvalues i=1/`counter' {
		mat graph[`i',1] = `bw_`i''
		mat graph[`i',2] = `b_`i''
		mat graph[`i',3] = `lci_`i''
		mat graph[`i',4] = `uci_`i''
		
	}

	mata : st_matrix("graph", sort(st_matrix("graph"), 1))
	preserve
	clear 
	svmat graph
	* graphs
	local line_2 = round(`x_line2',0.001)
	*attacks
	graph tw (sc graph2 graph1) (rcap graph4 graph3 graph1, lcolor(dkgreen)), title("") subtitle("") ///
	xtitle("Margin Victory") ytitle("Point Estimate") xli(`bw') xlabel(#8) xscale(range(`bw_half' `bw_double')) legend(off)  graphregion(fcolor(white))  yline(0)
	
	*Save graph
	graph export "Fig_E1_`y'.pdf", as(pdf) replace
	
	restore
	}
	}	
	

	
	cd "$dir"
	use "Data\cand_level_persist_rep.dta",clear
	cd "$dir\Figures"
	**Here is where the running variable goes.
	foreach x in margin_victory {
	**Here is where the outcome goes.
	*foreach y in $outcomes {
		foreach y in $outcomes5 {
		************
		* Linear
		*************
		*mserd approach, linear (get optimal bw)
		rdrobust `y' `x', all vce(cluster muni_code) p(1) level(95) 

		local bw=round(`e(h_r)',0.001)
		local bw_double = `bw'*2
		local bw_half = `bw'/2
		local rho =  `e(h_r)'/`e(b_r)'
		local counter=1
		local b_`counter' = round(`e(tau_cl)',0.001)
		local uci_`counter' = round(`e(ci_r_rb)',0.001)
		local lci_`counter' = round(`e(ci_l_rb)',0.001)
		local bw_`counter' = round(`e(h_r)',0.001)


		*According to Cattaneo, Idrobo and Titiunik, 
		*we should do double of the optimal bandwidht.
		local j = 1
		while `j'>0 {
			local counter = `counter'+1
			
			*CTT approach
			local b = `bw_double'/`rho'
			rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')
			
			local b_`counter' = round(`e(tau_cl)',0.001)
			local uci_`counter' = round(`e(ci_r_rb)',0.001)
			local lci_`counter' = round(`e(ci_l_rb)',0.001)
			local bw_`counter' = round(`e(h_r)',0.001)
			local bw_double = `bw_double'- (0.005)
			
			if `bw_double'<=`bw_half' {
				local j = 0
				
				* CTT approach
				local counter = `counter'+1
				local b = `bw_double'/`rho'
				rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')

				local b_`counter' = round(`e(tau_cl)',0.001)
				local uci_`counter' = round(`e(ci_r_rb)',0.001)
				local lci_`counter' = round(`e(ci_l_rb)',0.001)
				local bw_`counter' = round(`e(h_r)',0.001)
					
			
			}
			else {
				local j = 1
			}   
	}


	matrix graph = J(`counter',4,.)
	forvalues i=1/`counter' {
		mat graph[`i',1] = `bw_`i''
		mat graph[`i',2] = `b_`i''
		mat graph[`i',3] = `lci_`i''
		mat graph[`i',4] = `uci_`i''
		
	}

	mata : st_matrix("graph", sort(st_matrix("graph"), 1))
	preserve
	clear 
	svmat graph
	* graphs
	local line_2 = round(`x_line2',0.001)
	*attacks
	graph tw (sc graph2 graph1) (rcap graph4 graph3 graph1, lcolor(dkgreen)), title("") subtitle("") ///
	xtitle("Margin Victory") ytitle("Point Estimate") xli(`bw') xlabel(#8) xscale(range(`bw_half' `bw_double')) legend(off)  graphregion(fcolor(white))  yline(0)
	
	*Save graph
	graph export "Fig_E3_`y'.pdf", as(pdf) replace
	
	restore
	}
	}	
	
	cd "$dir"
	use "Data\cand_level_persist_rep.dta",clear
	cd "$dir\Figures"
	**Here is where the running variable goes.
	foreach x in margin_victory {
	**Here is where the outcome goes.
	*foreach y in $outcomes {
		foreach y in $outcomes4 {
		************
		* Linear
		*************
		*mserd approach, linear (get optimal bw)
		rdrobust `y' `x', all vce(cluster muni_code) p(1) level(95) 

		local bw=round(`e(h_r)',0.001)
		local bw_double = `bw'*2
		local bw_half = `bw'/2
		local rho =  `e(h_r)'/`e(b_r)'
		local counter=1
		local b_`counter' = round(`e(tau_cl)',0.001)
		local uci_`counter' = round(`e(ci_r_rb)',0.001)
		local lci_`counter' = round(`e(ci_l_rb)',0.001)
		local bw_`counter' = round(`e(h_r)',0.001)


		*According to Cattaneo, Idrobo and Titiunik, 
		*we should do double of the optimal bandwidht.
		local j = 1
		while `j'>0 {
			local counter = `counter'+1
			
			*CTT approach
			local b = `bw_double'/`rho'
			rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')
			
			local b_`counter' = round(`e(tau_cl)',0.001)
			local uci_`counter' = round(`e(ci_r_rb)',0.001)
			local lci_`counter' = round(`e(ci_l_rb)',0.001)
			local bw_`counter' = round(`e(h_r)',0.001)
			local bw_double = `bw_double'- (0.005)
			
			if `bw_double'<=`bw_half' {
				local j = 0
				
				* CTT approach
				local counter = `counter'+1
				local b = `bw_double'/`rho'
				rdrobust `y' `x', all vce(cluster muni_code) level(95) h(`bw_double') b(`b')

				local b_`counter' = round(`e(tau_cl)',0.001)
				local uci_`counter' = round(`e(ci_r_rb)',0.001)
				local lci_`counter' = round(`e(ci_l_rb)',0.001)
				local bw_`counter' = round(`e(h_r)',0.001)
					
			
			}
			else {
				local j = 1
			}   
	}


	matrix graph = J(`counter',4,.)
	forvalues i=1/`counter' {
		mat graph[`i',1] = `bw_`i''
		mat graph[`i',2] = `b_`i''
		mat graph[`i',3] = `lci_`i''
		mat graph[`i',4] = `uci_`i''
		
	}

	mata : st_matrix("graph", sort(st_matrix("graph"), 1))
	preserve
	clear 
	svmat graph
	* graphs
	local line_2 = round(`x_line2',0.001)
	*attacks
	graph tw (sc graph2 graph1) (rcap graph4 graph3 graph1, lcolor(dkgreen)), title("") subtitle("") ///
	xtitle("Margin Victory") ytitle("Point Estimate") xli(`bw') xlabel(#8) xscale(range(`bw_half' `bw_double')) legend(off)  graphregion(fcolor(white))  yline(0)
	
	*Save graph
	graph export "Fig_E3_`y'.pdf", as(pdf) replace
	
	restore
	}
	}